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Abstract 

Collider, direct and indirect searches for dark matter have typically little or no sensitivity to 
weakly interacting massive particles (WIMPs) with masses above a few TeV. This rather unex¬ 
plored regime can however be probed through the search for distinctive gamma-ray spectral features 
produced by the annihilation of WIMPs at very high energies. Here we present a dedicated search 
for gamma-ray boxes - sharp spectral features that cannot be mimicked by astrophysical sources 
- with the upcoming Cherenkov Telescope Array (CTA). Using realistic projections for the instru¬ 
ment performance and detailed background modelling, a profile likelihood analysis is implemented 
to derive the expected upper limits and sensitivity reach after 100 h of observations towards a 2° x 2° 
region around the Galactic centre. Our results show that CTA will be able to probe gamma-ray 
boxes down to annihilation cross sections of 1CV 27 — 10 -26 cm 3 /s up to tens of TeV. We also identify 
a number of concrete particle physics models providing thermal dark matter candidates that can 
be used as target benchmarks in future search campaigns. This constitutes a golden opportunity 
for CTA to either discover or rule out multi-TeV thermal dark matter in a corner of parameter 
space where all other experimental efforts are basically insensitive. 
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I. INTRODUCTION 


Among the myriad of dark matter candidates proposed, weakly interacting massive parti¬ 
cles (WIMPs) still stand among the best motivated possibilities [IHJj. Under very plausible 
conditions, the freeze-out at very early times of the WIMP number density EHZI leads to 
a relic abundance which is in qualitative agreement with the observed dark matter abun¬ 
dance, Odm^ 2 = 0.1198 ± 0.0015 [8], provided the thermally averaged annihilation cross 
section times the relative velocity of the dark matter particles at the time of freeze out is 
(av) ~ 1 pb • c. The dark matter mass, however, cannot be predicted in this framework. 

Most WIMP searches concentrate in the mass window between a few tens and a few 
hundreds of GeV, motivated by the possible relation of the dark matter problem with the 
stabilization of the electroweak scale against quadratically divergent quantum corrections. 
This is the case, for example, of the lightest neutralino in supersymmetric scenarios with 
R-parity conservation [;9j or the lightest Kaluza-Klein particle in scenarios with universal 
large extra dimensions ra- However, on model independent grounds, the theoretically 
allowed window for WIMP masses is much wider, and could be as large as ~ 118 (84) TeV for 
Majorana (Dirac) particles from the requirement of partial wave unitarity in the annihilation 
process [lTj (using the latest measurement of the relic abundance [8]), or as low as ~ 2 GeV, 
the renowned Lee-Weinberg limit for weakly interacting fermions |12j . or even lower for 
scalar dark matter [TB]. A complete investigation of the WIMP signals should therefore also 
include the high- and low-mass regions of the theoretically allowed parameter space. 

In this work we focus on the search for heavy dark matter particles. This mass regime is 
especially challenging for the Large Hadron Collider (LHC), which will not be able to discover 
signals of WIMPs heavier than about 1.5 TeV after running at 14 TeV centre-of-mass energy 
with 300 fb -1 of collected data mm • Future machines, such as a high-energy upgrade 
of the LHC with 33 TeV centre-of-mass energy [TB] or a 100 TeV proton-proton collider 
will be able to extend the reach to higher masses, but not any time soon. The high-mass 
regime is also challenging for direct detection experiments, which reach their maximum 
sensitivity for WIMP masses of several tens of GeV. For instance, the sensitivity of LUX 
|18J (XENON100 [19] ) is best at masses of ~ 33 (55) GeV under typical assumptions, and 
degrades quickly towards high masses. Nevertheless, the impressive limits on the scattering 
WIMP-nucleon cross section set by current experiments probe very heavy WIMPs in some 
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concrete scenarios (see e.g. 12®. 

Indirect search experiments, and more specifically searches using gamma rays and neutri¬ 
nos as messengers, may provide a very powerful probe for heavy dark matter. Some gamma- 
ray and neutrino telescopes are designed for the study of high-energy astrophysical phenom¬ 
ena and have their optimal sensitivity in the multi-TeV energy range, thus making them 
unique instruments to search for signals from the annihilation of heavy dark matter particles. 
There are several experiments currently operating which are probing WIMP masses above a 
few TeV, such as the gamma-ray instruments Fermi-LAT [2T] [22], H.E.S.S. [23, 2Tj, MAGIC 
[25] and VERITAS [2Bj, as well as the neutrino telescopes Super-Kamiokande 123 EH], Ice- 
Cube [2SSQ], ANTARES [31] and Baksan [32]. Unfortunately, the current limits are, for 
most scenarios, still far away from the required sensitivity to observe signals from heavy 
WIMP annihilations. 

A notable exception is provided by scenarios generating gamma-ray boxes [33]. This sharp 
feature in the gamma-ray energy spectrum cannot be mimicked by any known astrophysical 
process, so the observation of a gamma-ray box would constitute an unambiguous signal 
for dark matter annihilations. Moreover, in contrast to annihilations generating gamma- 
ray lines [33R36] or internal bremsstrahlung [371 [38], which have cross sections suppressed 
with respect to the thermal value by a factor 0(a 2 ) and 0(a), respectively, with a the 
electromagnetic coupling constant, the gamma-ray box has an a priori unsuppressed cross 
section and can therefore generate an intense, and at the same time very characteristic, signal 
in the gamma-ray spectrum. In fact, in the most favourable scenarios, current instruments 
are sensitive to gamma-ray boxes generated by annihilations of WIMPs with masses as large 
as ~ 10 TeV [39]. Motivated by the excellent sensitivity of current instruments to gamma- 
ray boxes, we investigate in this paper the prospects to observe this sharp feature with the 
future Cherenkov Telescope Array (CTA) [TOl HI]. 

The paper is organised as follows. In Section [TT] we review the most salient aspects of 
gamma-ray boxes and introduce a series of benchmark scenarios for our analysis. In Section 


i Id we describe our methodology to calculate the CTA reach and in Section IV we present 
our results. Lastly, in Section [V] we present our conclusions. 
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II. GAMMA-RAY BOXES FROM DARK MATTER 


The search for gamma-ray boxes originated in dark matter annihilations or decays has 
been proposed and discussed at length elsewhere [33, 33]. Here, we simply recall the essential 
phenomenology to compute the gamma-ray flux at Earth and outline concrete dark matter 
models for illustration purposes. 


A. Phenomenology 


Gamma-ray boxes are the result of dark matter cascade annihilations or decays into 
(pseudo-)scalar intermediate particles that are short-lived and decay into photons. Consider 
for concreteness the case of Dirac dark matter particles x annihilating into two scalars (j) 
which decay into a pair of photons each: XX — )• 00 — >• 4y (for variations on this scenario 
we refer the reader to Refs. [33] [33]). The scalars are sufficiently short-lived such that their 
decays occur basically in the same position as the dark matter annihilations - in practice, 
since we are interested in regions of the Galaxy spanning hundreds of parsec, decay times 
below about 10 2 — 10 3 yr are sufficiently short for our purposes. In the rest frame of the 
parent scalar, the photons are monochromatic with energy E' = m^j 2. After boosting 
to the lab frame, their energy E 7 depends on the emission angle and lies in between the 
kinematic limits 



where the degeneracy parameter 5ij = 1 — m 2 /m 2 has been used and E + (EL) corresponds 
to a photon emitted in the forward (backward) direction. Since the intermediate state is 
a scalar, the angular distribution of the photons is isotropic (in the rest frame of 0) and 
all lab-frame energies between E + and EL are equally populated, giving rise to a photon 
spectrum with a characteristic box shape: 


iwr^E e{E - E - )e{E+ - E) - (2) 

In the formula above, A E = E + — E_ is the width of the box and Q(x) is the Heaviside 
step function. The photon spectrum at Earth induced by such injection spectrum reads 


d$ dm _ ( <Jv)o 7 diV 7 
dE 7 I6irm^ dE 7 


with 


J :itil 




ds Pdm 


(3) 
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for symmetric Dirac dark matter particles, where (av) q 7 is the present-day thermally av¬ 
eraged annihilation cross section of the cascade \X —> </></> —» 4y, J ann represents the usual 
J-factor, Afi is the observed held of view, s parameterizes the direction down the line of 
sight defined by Galactic coordinates (£, b ) and pdm is the dark matter density. In this work, 
we shall focus on a 2 ° x 2 ° held of view around the Galactic centre (i.e. \l\ < 1 ° and | 6 | < 1 °) 
and consider different dark matter prohles: Einasto [42, 33] with scale radius r s = 20kpc and 
shape parameter a = 0.17, Navarro-Frenk-White (NFW) [441 145 ] with r s = 20kpc, a cuspy 
generalised (a, (3, 7 ) prohle [441146] with r s = 20 kpc and (a, (3, 7 ) = (1,3, 1 . 2 ) (motivated 
by the gamma-ray excess towards the Galactic centre m ) and a cored generalised (a, (3, 7 ) 
prohle with r s ~ 4.4 kpc and (ct,/3, 7 ) ~ (2.9,2.5,0) (motivated by recent high-resolution 
simulations including baryons [45]; we use halo B specihed in Tab. 3 in that reference). All 
prohles are normalised to p 0 = 0.4GeV/cm 3 at the position of the Sun, adopting a distance 
to the Galactic centre Ro = 8.5 kpc. For reference, the J-factors towards our 2° x 2° target 
region read J ann = 1.02 x 10 22 ,6.34 x 10 21 ,3.63 x 10 22 ,1.60 x 10 20 GeV 2 sr/cm 5 for Einasto, 
NFW, cuspy and cored prohles, respectively. 

The hux in Eq. (J3| assumes no attenuation of the photons in their route to Earth. 
However, as pointed out right after the discovery of the cosmic microwave background (CMB) 
[49H5T] , the universe is not transparent to high-energy gamma rays, their mean free path due 
to pair production on background photons reaching a minimum of ~ 8 kpc for E 7 ~ 2000 TeV 
(e.g. H). Since our region of interest is the centre of the Milky Way, at around 8 kpc 
from us, and we are in fact interested in energies up to hundreds of TeV, it is pertinent to 
calculate the optical depth of gamma rays due to pair production on the interstellar radiation 
held. The latter is composed of starlight (mostly ultraviolet), dust-reprocessed emission 
(mainly infrared) and the CMB. The attenuation due to starlight and dust emission becomes 
relevant at photon energies of a few tens of TeV, while the influence of CMB kicks in at 
around 200 TeV. We compute the optical depth r due to pair production following standard 
techniques (see e.g. [53]) and using the interstellar radiation held as found in [541 55] at 
Galactocentric radii R = 0,4, 12 and 16 kpc. Taking account of photon attenuation makes 
the J-factor depend on energy: 

JLn(Ey)= [ dSl [ dse-^'^yL, (4) 

J A J los 

where the optical depth t(Ej,s,Q) is computed along the line of sight. For our 2° x 2° 
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region of interest and the Einasto profile, we obtain JJ ni1 / J ann ~ 1, 0.99, 0.75, 0.35 for E 1 = 
10 GeV, 1 TeV, 100 TeV, 1000 TeV, respectively. 

B. Concrete models 

Physical models that feature boxed-shaped spectra are usually based on the breaking of 
global symmetries and involve Dirac dark matter [39J lobl iGTT] (see also [6Tj for a model with 
chiral fermion dark matter). We first discuss the setup developed in Refs. [39l 1561157] . where 
the Peccei-Quinn mechanism [62] is invoked and the resulting pseudoscalars play the role of 
the intermediate states (j) referred to in the previous subsection. The Lagrangian £ of this 
model reads 

£ = £sM + + Ant j (5) 

where £sm is the Standard Model Lagrangian, which includes the potential for the Higgs 
held, whereas 


, ( 6 ) 

£s = d'tS&S* - m 2 s \S\ 2 - A s |5 | 4 + Qm' s 2 S 2 + c.c.) , (7) 

4„t = -A x ( SxPlX + S’xPrX ) - 2X h , s \H\ 2 \S\ 2 ( 8 ) 

are, respectively, the parts of the Lagrangian involving only the held y, only the held S 
and the interaction Lagrangian. We further assume that the held S acquires a vacuum 
expectation value, so that the held S can be decomposed as 5 = (S) + ^(s + ia). The 
interaction Lagrangian in Eq. (J 8 ]) then generates a mass term for the Dirac held y, as well 
as a mixing term between the CP-even scalar s and the Standard Model Higgs boson. This 
mixing term, in particular, permits the thermalisation of all new states with the Standard 
Model plasma in the early universe. 

If the mass term m'g vanishes, the Lagrangian displays a global U( 1)pq symmetry under 
the held transformations y —> e* 75Q y and S — y e 2ia S , while all the Standard Model helds 
remain invariant. This symmetry is broken in the vacuum, however a residual Z 2 symmetry 
remains, under which y is odd and all other helds are even. The exact Z 2 symmetry then 
guarantees the stability of y, which constitutes our dark matter candidate. The observed 
dark matter abundance, moreover, can be explained with the freeze-out of y from the thermal 
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plasma due to the annihilations XX aa, ss i sa i when kinematically allowed. We note that 
the held a, being the Goldstone boson resulting from the spontaneous breaking of the global 
{7(1)pq symmetry, remains massless and is therefore also absolutely stable, thus constituting 
a candidate for dark radiation |63j |6T |. 

Here, we instead allow the mass term m'g to be non-vanishing, thus breaking explicitly 
the U( 1)pq symmetry. This variation does not affect the stability of the dark matter particle, 
however it introduces a mass for the pseudoscalar, m a = m' S) which can then decay into 
Standard Model particles. In particular, the pseudoscalar could decay into two electroweak 
gauge bosons if the model is extended with new heavy fermions with electroweak charges 
that couple to the complex scalar S, thus generating the following term in the interaction 
Lagrangian via the anomaly loop: 


ini ^ ^ 
i= 1,2 


CjCCj 

87TZU 


aFl.F 


l [lls 


[IV 


(9) 


and which induces the decays a —>■ 77, 7 Z, ZZ, W + W~. Here, <7 are constant parameters 
which depend on the assigned charges of heavy fermions in the anomaly loop, the parameters 
oti are given by cq = gf /An where g 12 are the coupling constants of U( 1)y and SU( 2) L , 
respectively, and finally F /JU = e^ vpa F pa /2 is the dual held strength tensor. 

The interactions described by the Lagrangian in Eq. ([ 8 ]) induce the dark matter annihila¬ 
tions XX aa j ss , sa. While s decays producing a continuum of gamma-rays, through the 
mixing with the Higgs boson, the pseudoscalar a can decay into one or two photons, hence 
producing a gamma-ray box. In order to assess the prospects to observe signals of this model 
we consider three different benchmark scenarios taken from Ref. |S9J, which assume m s = m a 
and that we label Al, A 2 and A3, corresponding to narrow, intermediate and wide boxes, 
respectively (see Tab. [Tj). I 11 all three cases the ratios for XX aa an d XX a A needed 
to produce a box spectrum, sum up to 50% or more. The decay of a into vector bosons 
depends on the constants cy we adopt ci/c 2 = 3 for which the branching ratio BR(a —» 77 ) 
is 100 % for m a < 2 mw and gradually decreases to 20 % at high pseudoscalar masses (for 
exact expressions see Ref. [32j). With these figures, the annihilation cross section in Eq. ([3]) 
reads 

(° v )V = ( av )o (bR(xx aa) + ^ BR(a 77 ), ( 10 ) 

where (crv) 0 is the present-day total annihilation cross section. The ratio of (crv) 0 to the 
corresponding value at thermal decoupling (cru) t h is shown in Tab. [Tj We compute the thermal 
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average in the non-relativistic limit as in Ref. [39]; we have however verified that a full 
treatment [63, US] leads to changes of less than 10% on the ratio (av) 0 /(av) t h- In order to 
test thermal candidates with the observed relic abundance, we shall take the canonical value 
(<7u) t h = 6 x 10 _ 26 cm 3 /s. The factor 2 with respect to the usual 3 x 10 _ 26 cm 3 /s is because 
the models considered here involve Dirac dark matter, not Majorana. 

Additionally, we analyse the framework presented in [o 8 l 16(11 . designed to interpret the 
positron excesses observed by PAMELA im Fermi-LAT [ 68 j and AMS-02 [69] while circum¬ 
venting the non-observation of large gamma-ray fluxes due to inverse Compton scattering 
[59] or a large antiproton flux. In these models, dark matter couples to a scalar s and 
a pseudoscalar a as in the previous framework, but now the s-wave annihilation channel 
XX y as is dominant [58]. Therefore, BR(xx — >• as) ~ 1 and (cru)o = ( crv)th as indi¬ 
cated in Tab. |T] for the benchmark scenario labelled B. Again, we take the canonical value 
(cry) t h = 6 x 10 -26 cn^s' 1 . Notice that this framework aims at explaining experimental data 
and therefore the masses of the involved particles are restricted to the ranges m x ~ (P(TeV), 
m s <C m x and m a = 360 — 800 MeV, giving rise to wide boxes only. For definiteness, we 
assume m a /m x = m s /m x = 0.1 throughout (note that the phenomenology of wide boxes is 
insensitive to the exact mass ratio, see [33] [39]). Furthermore, the decay of a is primarily 
leptonic leading to a branching ratio BR(a —> 77 ) ~ 10 -3 , which is however compensated 
by a Sonnnerfeld enhancement 1 B$ ~ 10 3 on the present-day annihilation rate [58]. In the 
end the annihilation cross section entering Eq. (|3]) is 

(<™)o 7 = B s (au) 0 BR( ^ 2 ^ aS ^ BR(a 77) . (11) 


The two classes of models presented above will serve as benchmarks for the following 
discussions, where Eqs. (10) and (cf. also Tab.[IJ will be used as target velocity-averaged 
cross sections to test the performance of the CTA. 


1 In principle, also scenarios A give rise to Sommerfeld enhancement, but to the best of our knowledge no 
previous study has addressed this behaviour. A full treatment lies outside of the scope of the current 
paper and is left for future work. 






m a /m x 

BR(xX -> as) 

BR(xX -> aa) 

BR(XX -> ss) 

(av) 0 /(av) th 

References 

Al (narrow) 

0.999 

0.99 

2 x 10 -3 

4 x 10 -3 

0.13 

E&HSHEZI 

A2 (intermediate) 

0.9 

0.64 

0.12 

0.24 

0.76 

[391 EH EH 

A3 (wide) 

0.1 

0.25 

0.25 

0.50 

0.96 

m ESI EH 

B (wide) 

0.1 

1 

0 

0 

1 

EM] 


TABLE I: Details of the concrete models based on Refs. [39 s 56i E] (Al, A2, A3) giving rise to 
narrow, intermediate and wide boxes and on Refs. [5&H60] (B) where only wide boxes are produced. 
In all cases we have assumed m s = m a - 


III. METHODOLOGY 

In order to derive the prospects for CTA to detect or rule out gamma-ray boxes we need 
a detailed account of the instrument properties, the expected background and a statistical 
procedure to convey our results. These three topics are now discussed in turn. 


A. CTA performance 


The upcoming Cherenkov Telescope Array (CTA) j30j 1H| will consist of two sites, one in 
each hemisphere, equipped with three different kinds of telescopes in order to cover a wide 
energy range from 10 GeV to about 100 TeV. CTA is expected to outperform existing arrays 
by one order of magnitude |41j . reaching differential sensitivities down to a fewmilliCrab 
at 1 TeV after 50 h of observations dSE]. Detailed information on Monte Carlo design 
studies can be found in Ref. DU- In our work we implement the performance of the southern 
candidate array I, which is a well-studied, balanced configuration consisting of 3 large, 18 
medium and 56 small size telescopes. Above a few TeV the effective area of this array exceeds 
10 6 m 2 and the energy resolution is better than 10%. We use throughout the detailed effective 
area A e g and energy resolution a e as found using the baseline MPIK analysis (cf. Figs. 15 
and 18 in Ref. mi). Note that the effective area in Ref. HU is given at the analysis level so 
that we take a photon acceptance e 7 = 1. Our study does not explore the angular resolution 
of the instrument, but instead focuses on a 2° x 2° region around the Galactic centre (as 


mentioned in Sec. IIA) for which an observation time At = 100h is reasonable given the 


high astrophysical interest of this region. We further assume an electron acceptance e e = 1 
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and a proton acceptance of e p = 0.01 — 0.20 (see Sec. IllB for a full expression). Finally, let 
us point out that recently an updated study of the expected performance of CTA has been 
made available (72j : we also implement this updated performance for which effective area 
and energy resolution are better than those of Ref. ra and the proton acceptance reads 


instead e p = 0.001 — 0.02 (again, see Sec. IllB for further details). The performance based 
on Ref. [[2] is used throughout the work and we comment on the results obtained with the 
updated performance [72] in Sec. IV (cf. in particular Figs. [2] and [5]). 


B. Backgrounds 

The search for dark matter with imaging air Cherenkov telescopes is usually hindered by 
three sorts of background: 

• First, hadronic cosmic rays (mainly protons) initiate showers whose Cherenkov radia¬ 
tion may be misread as coming from a photon-initiated shower. It is actually possible 
to discriminate hadronic and electromagnetic cascades - typical hadron acceptances 
vary between 1% and 20% for a 70% gamma-ray acceptance [73]. However, cosmic 
rays are still a sizeable background since they largely outnumber photons. 

• Second, showers induced by electrons or positrons are indistinguishable from those 
triggered by photons at the same energy (apart from a slightly later shower maximum 
on average for photons ES) and constitute an important background as well. 

• Third, photons from astrophysical or other origin ultimately limit the sensitivity of 
dark matter searches and can only be mitigated with a wise choice of the target field 
of view and analysis technique. 

In the following we describe in detail the parameterizations used in our analysis for each 
type of background. Let us start with cosmic-ray protons. We follow Ref. E5j where it was 
found that the power law 

= 8.73 x 10~ 6 ( TeV -1 cm -2 s -1 sr -1 (12) 

dFpdfl V TeV / 

provides a good fit to a wide array of data prior to 2002 from around 10 GeV up to knee 
energies. Including the latest proton data ca changes only slightly the best-fit power 
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law, so we use the above parameterization throughout. Notice that for primary energies 
E 0 > 100 GeV the Cherenkov yield of a photon shower is a factor p ~ 2 — 3 larger than that 
of a proton shower of the same energy, and both yields are roughly proportional to E 0 |77j. 
Therefore, a shower triggered by a proton with energy E p has the same yield as a photon- 
initiated shower of a smaller primary energy IT, ~ E p / p. Actually, there are several shower 
parameters other than the Cherenkov yield used to discriminate photons from hadrons m, 
but, in the absence of detailed information about the typical energies of misreconstructed 
protons, we use the relation above with p = 3 in computing the background flux due to 
protons. 

The all-electron spectrum has also been precisely determined, namely by AMS-02 na 
below several hundred GeV and by H.E. s.s. m at TeV energies. AMS-02 data follows 
essentially a power law 


d 2 ^ 

dE e dil 


9.93 x 10~ 9 


/ 77 \ 

f tjt-A j TeV _1 cm~ 2 s _1 sr _1 


(13) 


for energies above 30.2 GeV, while H.E.S.S. data yields 


d 2< f>g e 

dE e dfl 


1.17 x 10“ 8 


77 \ —3.9 

\ rp u -1 -2 -1 -1 

——- lev cm s sr 

TeV J 


(14) 


In order to ensure a smooth transition between the two regimes the total all-electron flux 
reads d 2 $ e /d£ e dfi = [(d 2 ^/dE e dQ)~ 2 + (d 2 ^ e /dE e dQ)- 2 ]~ 1 / 2 . 

Finally, the background due to gamma rays of astrophysical origin depends heavily on the 
chosen target held of view. In our 2° x 2° region of interest there are two main gamma-ray 
components at high energies. On the one hand, the Galactic ridge emission spanning the 
region |£| < 0.8° and |6| < 0.3° has been detected by H.E.S.S. [79] and is well described by 
the power law 

rl 2 d) / 77 \ —2.29 

d^ = L 73 xl 0 "(sv) ,;v (15) 

On the other hand, data from Fermi-LAT and H.E.S.S. show a point source coincident with 
the Galactic centre (to within the resolution of the instruments) whose energy spectrum in 
the range 5 — 100 GeV is well fitted by [80] 

d<h le f E \" 2 - 68 

—^ = 1.11 x 10- 12 TeV- 1 cm- 2 s- 1 (16) 

d E 1 \TeV / v ; 
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and in the range 160 GeV — 30 TeV by |81 


(Mg" 19 ( e~ 

- 7,s = 2.34 x 10~ 12 —A- 
d E~ \ TeV 


-2.25 


-l„„-2„-l 


TeV cm 


(17) 


A smooth transition between low and high energies is achieved using d < h 7igc /di? 7 = 

[(dt^/dB,) 5 + (d^/dE,) 5 ]'/'. 

The overall background flux then reads 

d$bk g s . „ ( , d 2 $„ , „ , d 2< h P , d 2< h^.,, r A d<f> 


dA, 


( d 2 <h d 2 <h 

(E 7 ) = Afi 2x2 ( P E 7 ) + e 


dEAQ 


{E~f) c 


-7.gr 

7 dE 7 dfi 


(E*,) )+e 


7>g c 


d E, 


(£7) , 


7 

(18) 

where Ah2 2x2 = 1.22 x 10~ 3 sr is the solid angle of our 2° x 2° region of interest. Notice that, 
in order to be as conservative as possible, we upscale the flux found in the Galactic ridge 
(\£\ < 0.8° and |6| < 0.3°) to match our target field of view assuming the same output per 
unit solid angle. Furthermore, we have checked that our cosmic-ray background rate matches 
the one reported in Ref. HU (see Fig. 16 therein) within a factor 1.5 — 3 if we take e p (E 7 ) = 
0.01 + 0.02(E 7 /20 TeV) 1 ' 4 , i.e. a benchmark 1% proton acceptance below ~ 1 TeV that 
grows up to around 20% at 100 TeV. This is roughly in line with Ref. [82], where a full CTA 
simulation is used to set the cosmic-ray background rate. We repeated the same procedure 
with the updated CTA performance [72], obtaining e p {E 1 ) = 0.001 + 0.002(E 7 /20TeV) 1 ' 4 . 

The breakdown of the background in Eq. (18) is shown in the left plot of Fig. [I] Clearly, 
below a few hundred GeV cosmic-ray electrons and positrons represent the largest back¬ 
ground, whereas at TeV energies and above the diffuse emission from the Galactic ridge 
and eventually protons overshadow all other components. The right panel of Fig. [l] shows 
a side-by-side comparison of the overall background and the signal from Eq. (J3]) for the 
benchmark values (av) q 7 = 10 _24 cm 3 /s, m x = 20 TeV and m^jm x = 0.999, 0.9, 0.1. 


C. Limits and sensitivity 

After outlining the properties of CTA and defining background and signal fluxes, we can 
now generate mock data for CTA. For definiteness, Nf, = 200 energy bins per decade are set 
up in the energy range E 1 ~ 40 GeV — 110 TeV. Assuming a generic flux d$ to t/d^ 7 , the 
expected number of counts n* xp in bin i of width A E t is simply 

»U = A t f d E f d E' R(E, E') A m (E') (V), (19) 

J t/ 7 
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FIG. 1: The background for CTA gamma-ray searches in a 2° x 2° region around the Galactic centre. 


The left panel shows the breakdown of the total background in Eq. (18) into protons, electrons and 
positrons, and gamma rays from both the Galactic centre source and the Galactic ridge region. The 
right panel compares the total background flux with narrow, intermediate and wide gamma-ray 
boxes induced by a dark matter scenario with an Einasto profile, (av) q 7 = 10 -24 cm 3 /s, m x = 
20TeV and m ( p/m x = 0.999, 0.9, 0.1, respectively. Both background and signal in the right panel 
are convoluted with the CTA energy resolution described in Sec. Ill A| Here we use the CTA 
performance of Ref. m- 


where R(E, E') is the energy resolution of the instrument, assumed Gaussian with standard 


deviation given by erg (cf. Sec. Ill A). The mock number of counts n Q bs is drawn from a 


Poissonian distribution with mean rii 


exp* 


Given a data set }, we make use of a profile likelihood analysis jHITI 157] in order 

to search for gamma-ray boxes with fixed mass configuration ( m x ,m </,). Our flux model 
dc&tot/d-E^ includes a power-law background (with normalisation A and spectral index —7) 
and a putative gamma-ray box (with strength set by S = (<7u)q 7 ), in a total of three 
parameters. For each bin i in the data set, the expected number of counts n* xp (A, 7, S) 


follows from Eq. (19), and the likelihood is written as a product of Poissonian probabilities: 

exp(-<„) 


£=n 


(n 


exp ) 


n 


( 20 ) 


obs' 


We use the sliding window technique [HIHEH] such that the product in the expression above 
is restricted to the bins inside a given energy window [E/ yA, Ey/e\. The quantity —2 In C is 
then minimised over (A, 7) in order to define the profile likelihood — 21n£ pro f as a function 
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of S, and over (A, 7, S) to find the global minimum (i.e. best fit) —2 In £bf attained at 
(Abf, 7bf, S'bf) • The best fit under the null hypothesis follows immediately from the profile 
likelihood as — 21n£° f = — 21n£ pro f(S' = 0). 

The aim of this paper is to infer the prospects for CTA in the search for gamma-ray 
boxes, namely upper limits in case of non-observation and sensitivity reach otherwise. For 
the upper limits, the starting point is to generate 300 sets of mock data assuming back¬ 


ground only, i.e. with d<£> t ot/dE 7 = d<f>bkg/dE 7 in Eq. (19). Next, for each data set and any 
given combination (m x , m^), we compute the profile likelihood and find the one-sided 95% 
confidence level (CL) upper limit S u \ > Sbf defined as —2 In £ prof (5 u i) = —21n£ b f + 2.71 
|90J. For the sensitivity reach, having chosen a specific mass configuration (m x ,m^), 
we start by drawing 300 realisations of mock data assuming signal plus background, 
i.e. d<Ftot/dE 7 = e 7 d<Fdm/dE 7 + d<Fbkg/dE 7 , where the signal is computed with a rela¬ 
tively small value of S (i.e. (< 7 u)q 7 ). We then determine the average of the test statistic 
TS = — 2 (In £[( f — In £ b f) over the 300 realisations and repeat the procedure with a higher 
value of S. The 5cr sensitivity reach S^ corresponds to the smallest annihilation cross 
section for which the average TS is larger than 23.7 (see e.g. [Si]). 

Let us finally comment on the energy windows used in our analysis. We set the central 
energy of all windows to E — E + (cf. Eq. ([Tj) ) and adopt three different widths e: 61.2 = 1.2, 
61.5 = 1.5 and e 2 = 2. These windows are chosen to test the robustness of our results and, at 
the same time, ensure that the background inside the window is well fitted by a power law. 
Specifically, we started by taking the 300 mock data sets generated assuming background 
only (see above), and fitted each data set to a pure power law using windows of constant 
widths (spanning the range e = 1.2 — 10) and centred across the full CTA energy range. 
For each window, we thus obtain 300 y 2 values for the best fit power laws and test whether 
this distribution is plausibly drawn from a y 2 distribution with the appropriate degrees of 
freedom (number of bins minus two fitted parameters). Windows for which the p-value 
of this distribution test is smaller than 0.01 are rejected. The window of width e 2 is in 
particular the largest window surviving this criterion, so we shall adopt it as our default 
window throughout the work, while still using ei . 2 and 61.5 for completeness. 
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FIG. 2: Expected CTA upper limits and sensitivity reach for gamma-ray boxes in a 2° x 2° region 
around the Galactic centre after 100 h of observations. The left, centre and right panels refer 
to narrow, intermediate and wide boxes, respectively. The violet solid lines and corresponding 
bands show the average and standard deviation of the 95% CL upper limits on the annihilation 
cross section over 300 mock data sets based on the CTA performance of Ref. m, while the 5 a 
sensitivities are indicated by the dashed lines. The dotted lines correspond to the average 95% CL 
upper limits obtained with the updated CTA performance m- Here we assumed an Einasto dark 
matter profile, a vanishing gamma-ray optical depth and a window width C 2 - The thick solid lines 
show the target cross sections corresponding to thermal relics in models Al-3 and B of Tab. [T} 


IV. RESULTS 

The typology of box spectra is rather rich, comprising line-like features at half the dark 
matter mass, flat shoulders up to the dark matter mass and every variant in-between these 
two extremes. As discussed in Sec. |TT[ the shape of the spectrum is solely defined by the mass 
ratio m ( p/m x . Here we focus on the three benchmark values m ( j ) /m x = 0.999, 0.9, 0.1 that 
give rise to narrow, intermediate and wide boxes (cf. right panel of Fig. 0- Rg.0 shows 
for each box case the 95% CL upper limits and 5 a sensitivity reach expected after 100 h 
of CTA observations towards a 2° x 2° region around the Galactic centre. For this figure 
we have used our baseline assumptions, namely the Einasto profile, a vanishing gamma-ray 
optical depth and the sliding window of width 62 - Generically speaking, the average upper 
limits cover annihilation cross sections of a few times 10 ~ 2 ' cm 3 /s or smaller for dark matter 
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FIG. 3: Expected CTA upper limits for narrow (left panel), intermediate (centre) and wide (right) 
boxes. In each panel we show the results for the three sliding window widths £ 1 . 2 , 61.5 and £ 2 . 
All unspecified parameters or assumptions are fixed to the ones adopted in Fig. [2] and the CTA 
performance of Ref. m is used. 


masses from about 100 GeV up to 10 — 20 TeV, while the expected statistical variance on 
these limits corresponds to a factor 2 — 3. The 5 a sensitivity reach instead is worse than the 
average 95% CL upper limits by a factor 3 — 5 across the full dark matter mass range. Note 
that the expected prospects for narrow, intermediate and wide boxes are very similar at face 
value. As noticed in Refs. [33, 39], this happens because, although narrow features are the 
most convenient to search for against a smooth background, wider spectra extend to higher 
energies where the background is dimmer, partially compensating for their relative softness. 
In Fig. [2j we also show the average 95% CL upper limits obtained with the updated CTA 
performance 1721. which are better than the limits based on Ref. HD by a factor rs-/ 8 (5.6) 
at low masses for narrow (wide) boxes and at most by a factor of ~ 1.6 above 4 TeV for all 
three box sizes. A full study of the implications of the updated CTA performance [72] is 
deferred for future work. 

We now proceed to check the robustness of the results against our baseline assumptions. 
The effect of the sliding window on the expected CTA upper limits is presented in Fig. [3j 
Clearly, a wise choice of the window will be crucial to make the best of the available gamma- 
ray data. This is important for intermediate and wide boxes, where the tightest window of 
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FIG. 4: Expected CTA upper limits for different dark matter halo profiles, and effect of including 
the optical depth, for narrow (left panel), intermediate (centre) and wide (right) boxes. In each case 
we show only the average of the 95% CL upper limits. All unspecified parameters or assumptions 
are fixed to the ones adopted in Fig. [ 2 ] and the CTA performance of Ref. mi is used. 


width ei . 2 fails to capture a significant part of the signal giving thus the weakest constraints. 
Such behaviour is also seen for narrow boxes at low energies since there the CTA energy 
resolution is poor and the expected signal is somewhat extended in energy. In all cases in 
Fig. [3j larger windows can improve the limits as much as a factor 25 provided the background 
is well described by a power law. We have further checked that our baseline upper limits 
do not change significantly for different choices of the number of bins N b , and they can be 
improved up to a factor 1.5 — 2 by shifting the window centre E. 

In Fig. [4] we show the average 95% CL upper limits obtained for different dark matter 
profiles. If the dark matter distribution in the inner Galaxy follows a cuspy profile with inner 
slope 7 = 1.2, then our baseline limits reported in Fig. [2] improve by roughly a factor 4, 
whereas a cored profile leads to a weakening of almost two orders of magnitude. Also shown 
in Fig. [4] is the effect of the gamma-ray optical depth which, according to the estimates at 


the end of Sec. II A, accounts for an absorption of 25% of the 100TeV photons from dark 
matter annihilation towards the Galactic centre and deteriorates the limits by around 20% 
at the highest masses considered. 

It is instructive to compare our CTA prospects with the current limits on gamma-ray 
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FIG. 5: Comparison of CTA expected upper limits for narrow boxes with current Fermi-LAT m 
and H.E.S.S. [23] line searches. The left panel shows the upper limits on the cross section into 
four photons (both Fermi-LAT and H.E.S.S. results were appropriately rescaled to be compared 
to narrow boxes, i.e. lines at half the dark matter mass with four photons per annihilation), while 
the right panel displays the flux upper limits for a monochromatic signal generated at a given 
energy. The Fermi-LAT limits correspond to 3.7 yr of data towards a circular 3° region around the 
Galactic centre with J ann = 1.39 x 10 23 GeV 2 sr/cm 5 for an NFW contracted profile, whereas the 
H.E.S.S. searches use 112 h observations towards a circular 1° region around the Galactic centre 
excluding |6| < 0.3° with J ann = 4.43 x 10 21 GeV 2 sr/cm 5 for an Einasto profile. For the CTA 
limits shown here we assumed the same observation time, profile, region of interest and J-factor as 
in H.E.S.S. [23] ? an d all other parameters or assumptions are fixed to the ones adopted in Fig. [2j 
As in Fig. [2j the violet solid lines and bands correspond to the average upper limits obtained with 
the CTA performance of Ref. |71] , while the dotted lines are for the updated performance [72]. All 
bounds in this figure are 95% CL upper limits. 

spectral features. Such comparison is done in the left panel of Fig. [5] using the latest line 
searches from Fermi-LAT [9T] and H.E.S.S. {2Tj, which have been rescaled from the usual 
monochromatic line at m x with two photons per annihilation to a narrow box centred at 
m x j 2 with four photons per annihilation. The plotted limits correspond to reasonably similar 
setups, namely a circular 3° region around the Galactic centre and an NFW contracted 
profile for 3.7yr of Fermi-LAT data, and a circular 1° region around the Galactic centre 
excluding \b\ < 0.3° and an Einasto profile for 112 h of H.E.S.S. observations. For the 
purposes of comparison, we have run a brand new set of CTA constraints matching the 


18 







observation time, profile, region of interest and J-factor of the H.E.S.S. limits. The resulting 
upper limits are overplotted in Fig. [5] for the window with width using both the CTA 
performance of Ref. CB and the updated one [723- Our analysis suggests that at energies 
below a few hundred GeV CTA will hardly supersede the current bounds provided by Fermi- 
LAT, whereas above 1 TeV it will improve upon H.E.S.S. by up to one order of magnitude. 

In principle, we could compare our results for narrow boxes to previous works deriving line 
prospects for CTA, e.g. [85, 92j. However, let us point out that these works precede the latest 
CTA design studies of Ref. j7T] (and [72]), used in our analysis. A direct comparison would 
therefore be somewhat misleading since the CTA performance (in particular, effective area) 
implemented in Refs. [85] [92] differs substantially from Ref. |7TJ (and [72]). Nevertheless, 
in order to be as transparent as possible and enable a straightforward comparison between 
different works in the literature, we provide in the right panel of Fig. [5] the flux upper limits 
for monochromatic features at given energy IT,, i.e. E x = m x for a line or IT, = m x j 2 for a 
narrow box. 

At this point we should comment on the validity and limitations of our results. All 
the upper limits shown are naturally based on a background-only hypothesis, so our cross 
section constraints scale as 1/yfKt. We have explicitly verified that this behaviour roughly 
holds true with an additional run adopting 50 h of CTA observation time. In general, the 
results presented here are certainly not the last word in terms of CTA prospects: they can 
be improved with an optimised morphological analysis as pursued e.g. in Refs. p32, 23], |MJ 
and they are limited by the (unknown) systematic level of the instrument as studied in 
Refs. [821 IS]- Following the same approach as in Refs. [82l IB], we find that the effect of 
systematics in our limits is practically negligible at high masses and sizeable at low masses, 
in agreement with those references. In any case, the prospects presented in Figs. [2] through 
[5] constitute the Erst step to probe gamma-ray boxes from dark matter annihilation over the 
coming decade. 

Finally, we turn to the discussion of the concrete dark matter models introduced in 


Sec. |II B| Let us stress that the use of these models simply serves as illustration for the 
reach of CTA and is not intended as a comprehensive model building analysis. In Fig. [2] we 
show the cross sections (crv) q 7 expected from the annihilation of the thermal candidates Al-3 
and B in Tab. |T[ For the models Al-3 there is a step-wise behaviour, because the branching 


ratio BR(a —>■ 77 ) falls sharply for m a > 2 rriw (cf. Sec. IIB) and because for fixed mass 
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FIG. 6: The ratio between the average 95% CL upper limit expected for CTA and the target 
cross section of a thermal dark matter candidate in models A1-3 and B of Tab. [TJ All unspecified 
parameters or assumptions are fixed to the ones adopted in Fig. [2] and the CTA performance of 
Ref. [7TJ is used. 


ratio m a is proportional to m x . In the case of models Al and A2, corresponding to narrow 
and intermediate boxes, respectively, the sensitivity of CTA is sufficient to probe thermal 
candidates for relatively low masses (below 1 TeV), a regime in which other instruments 
(e.g. Fermi-LAT) provide stronger constraints. The most promising case involves wide boxes 
as in models A3 and B, for which thermal relics can be tested up to several TeV. We recall 
that all the above is valid for the assumed (smooth) dark matter profile and concrete particle 
physics model; therefore for convenience of model builders we provide in Fig. [6] the results 
of Fig. precast in terms of the ratio of the average 95% CL upper limit to the cross section 
expected from a thermal relic. If we have underestimated the expected signal flux (due to 
the dumpiness of the dark matter distribution |95[96], Sonnnerfeld enhancement [97H99] or 
other) by a modest factor ~ 10, then the models Al-3 and B can be probed with CTA up 
to tens of TeV. 
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V. CONCLUSION 


Gamma-ray spectral features at high energies provide perhaps the best way to eventually 
detect dark matter in indirect searches. Until recently, most efforts in this respect focussed 
on the search for monochromatic lines and internal bremsstrahlung, both proposed long ago. 
Gamma-ray boxes constitute a third type of spectral features which opens the possibility to 
test a whole new class of thermal dark matter candidates. This paper is a first attempt at 
exploring the capability of future Cherenkov instruments such as CTA to detect or rule out 
gamma-ray boxes in the high-energy regime. We find that CTA has indeed fairly good dis¬ 
covery prospects for thermal relics in the multi-TeV range, especially so for models featuring 
wide boxes. In case no signal is found, a large portion of thermal dark matter models will be 
robustly ruled out up to tens of TeV. This is a golden opportunity for CTA since virtually 
no other instrument or technique can access this mass regime. The mass range above a few 
tens of TeV and up to the roof given by unitarity (mj m A 84 TeV [8], [llj for Dirac dark 
matter) is very challenging to probe even with gamma-ray features and will probably remain 
open for a very long time. 
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Erratum: On the sensitivity of CTA to gamma-ray boxes from multi-TeV 
dark matter 

Abstract Following the erratum of Ref. [3TJj, we update here the target cross sections 
corresponding to the particle physics models A1-3 presented in the published manuscript. 

The corrections presented in the erratum of Ref. [39j have an impact in the cross sections 
and branching ratios of the three benchmark models A1-3 discussed in Sec. 2.2. This update 
affects Tab. [T] and Figs. [2] and [6} In Tab. |TJ after the corrections, all models Al-3 feature 
BR(xx as ) — I; while the cross section ratio (crv)o/{av ) t h is 0.24, 1.02 and 1.01 for 
models Al, A2 and A3, respectively. Overall, this translates into a target cross section 
(cru)o 7 (cf. Eq. (2.10)) larger by a factor of 1.82, 1.51 and 1.39 for models Al, A2 and A3, 
respectively. Correspondingly, the thick solid lines labelled Al-3 in Fig. [2] are to be rescaled 
upwards by these factors and the lines labelled Al-3 in Fig. [6] are to be rescaled downwards 
by the same factors. In the process of spotting this issue, we have also found out that in the 
published version of the manuscript the plots of Fig. [2] are missing the benchmark model 
cross sections (thick solid lines). Therefore, for completeness, we provide in Figs. pFj and [ti 7 ] 
the updated and corrected version of Figs. [2] an di These changes make the CTA prospects 
to constrain models Al-3 slightly more optimistic than presented in the published version 
of the work. None of our conclusions changes. 
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FIG. 2': Expected CTA upper limits and sensitivity reach for gamma-ray boxes in a 2° x 2° region 
around the Galactic centre after 100 h of observations. For a full description, please refer to the 
caption of Fig. 2 in the published version of the manuscript. Note that here the benchmark cross 
sections corresponding to models Al-3 have been updated according to the corrections described 


in the current erratum. 
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FIG. 6': The ratio between the average 95% CL upper limit expected for CTA and the target cross 
section of a thermal dark matter candidate in models Al-3 and B of Tab. 1. For a full description, 
please refer to the caption of Fig. 6 in the published version of the manuscript. Note that here 
the benchmark cross sections corresponding to models Al-3 have been updated according to the 
corrections described in the current erratum. 
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